sommerfeld: fix subdivision cap, interpolation cache, and efld radius - #130

Description

Over a GN 2 Sommerfeld ground, nec2++ reports feed impedances that are
wrong by tens to hundreds of ohms, and on low antennas the feed resistance
comes out negative — a passive structure appearing to generate power. A
40 m-band dipole (37-3e) reports R = -60.3 Ω where the NEC-2 Fortran reports
+130.194 Ω.

Three independent defects in the Sommerfeld path produce this. Each is
corrected in its own commit and each was measured in isolation with the other
two applied, because their effects are not additive: the thin-wire radius
correction moves 37-3e by 2.9 Ω while the interpolation cache is still
defective, and by 0.03 Ω once it is fixed.

Type of Change

• Bug fix (non-breaking change that fixes an issue)

C-lineage scope

One of the three defects is not specific to nec2++. The NEC-2 Fortran
declares the Romberg subdivision limit as a routine-local constant:

      DATA NM,NTS,NX,N/65536,4,1,9/,RX/1.D-4/

nec2dx.f:7976 and nec2dxs.f:9068 both carry it. The original C translation
substituted the structure's segment-plus-patch count in its place, and every
C-lineage port inherited that substitution:

The same fix has been applied to nec2c and xnec2c; only nec2++ reaches
the bound on published models, because the other two defects here drive it
there. The remaining two defects are nec2++-only.

Implementation Details

src/nec_context.cpp — rom2 subdivision cap (b765112)

Both rom2 copies read m_geometry->n_plus_m() at the two sites the reference
reads NM: the minimum-interval epsilon and the step-halving bound. The
accuracy ceiling of the ground integration therefore scaled with model size
rather than with the sharpness of the integrand, so the smallest models over
ground received the coarsest integration — exactly where per-segment ground
interaction is strongest.

When the 5-point Romberg test fails at the finest permitted step, the warning
path falls through to loop10, which accumulates t20 — the extrapolation
that just failed its own error test. The corrected run terminates on the error
test RX = 1e-4 instead of on an exhausted counter.

Both the live copy under #if 1 and the disabled twin under #else are
patched in one commit, so the toggle cannot reintroduce a divergent bound.

src/c_ggrid.cpp — interpolation cache guard (d0ee6a6)

c_ggrid::interpolate caches the bivariate coefficients of the 4-by-4 grid
region around the previous request. Its opening guard had the two branches
exchanged relative to SUBROUTINE INTRP: a point at or above the region
origin recomputed unconditionally, while a point below the origin — the case
the reference forces to rebuild, since it carries no valid index — reused
whatever coefficients the previous call left behind whenever the stale indices
happened to satisfy the proximity test.

The routine therefore returned different values for the same (x, y)
depending on the order points were requested. Restored to the reference form:
a jump flag for the below-origin case, indices computed only in the other
branch, and one disjunction combining jump with the four index tests.

src/nec_context.cpp — thin-wire displacement radius (a89a6f9)

efld displaces the observation point off the wire axis by the wire radius,
dividing by the squared radial distance. The Sommerfeld section computed
norm() where math_util.h offers norm2() for the squared sum, and the next
statement takes a square root of it — so the displacement was ai/sqrt(|rho|)
instead of ai/|rho|. The comment on the defective line stated the correct
expression.

The overloads invite this: for three scalars norm() is the square root of the
sum of squares, while for a complex argument norm() is the squared modulus.
The other twelve call sites in the file were checked against that distinction
and each is correct as written.

Measured proof

Feed resistance in ohms, nec2++ -w, each defect toggled alone with the other
two corrected. Reference is nec2dxs, which nec2dx reproduces exactly:

Order dependence of the cache defect, isolated from every other difference by
driving nec2++ interpolation over a grid nec2c generated. 20301 points
swept ascending, then the identical points descending:

Across the 16 GN 2 decks of a five-implementation validation corpus
(nec2c, xnec2c, nec2++, nec2dx, nec2dxs), failing comparison points
on nec2++ edges fall from 489 to 208 with all three fixes applied.

Reproduction

A one-wavelength Beverage isolates the subdivision bound from every other
property of a model: eleven segments give a segment-plus-patch count of 11,
while the integrand at 0.0015 wavelengths height needs ns = 32, so the
substituted bound blocks every doubling past 16.

Before: 22 ROM2 -- STEP SIZE LIMITED lines on stderr; feed R = 604.74 Ω.
After: no warnings; feed R = 604.87 Ω against a reference of 604.924 Ω.

CM Low Beverage receiving antenna: 170 m of #14 wire 0.25 m over average
CM earth at 1.8 MHz, fed at one end and terminated at the far end into
CM 450 ohms. One wavelength long, eleven segments of 0.093 wl each, the
CM wire at 0.0015 wl.
CE
GW 1 11 0 0 0.25 170 0 0.25 0.0008128
GE 1
GN 2 0 0 0 13.0 0.005
LD 4 1 11 11 450 0
FR 0 1 0 0 1.8 0
EX 0 1 1 0 1 0
RP 0 10 37 1000 0 0 10 10
EN

The same antenna raised to 10 m settles at ns = 2 and its output is
byte-identical before and after the subdivision fix, so height alone drives
the difference. It is not inert to the other two defects: with the cache guard
reverted it moves 48 Ω in reactance, which is what makes it a control rather
than a dead deck.

Before: GW 1 11 0 0 10 170 0 10 0.0008128, feed 3250.9 + j615.17 Ω with the
cache guard defective.
After: 3264.5 + j567.17 Ω, reference 3262.85 + j567.229 Ω.

Testing

• Cross-engine validation against nec2dx and nec2dxs on 16 GN 2 decks
• Feed-point agreement to printed precision on 37-3e and 11-3a
• Order-dependence sweep over 20301 interpolation points, forward and
  reverse, on a grid generated by an independent implementation
• Upstream nec2++ test suite not run in this environment

The 37-3e, 11-3a and 20-40m_vert_sommerfeld_cliff decks are published
models; the Beverage pair was authored for this defect and is the minimal
witness.